Tuning the structural and dynamical properties of a dipolar Bose-Einstein 
condensate: Ripples and instability islands 
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It is now well established that the stability of aligned dipolar Bose gases can be tuned by varying 
the aspect ratio of the external harmonic confinement. This paper extends this idea and demon- 
strates that a Gaussian barrier along the strong confinement direction can be employed to tune both 
the structural properties and the dynamical stability of an oblate dipolar Bose gas aligned along 
the strong confinement direction. In particular, our theoretical mean-field analysis predicts the ex- 
istence of instability islands immersed in otherwise stable regions of the phase diagram. Dynamical 
studies indicate that these instability islands, which can be probed experimentally with present-day 
technology, are associated with the going soft of a Bogoliubov-de Gennes excitation frequency with 
radial breathing mode character. Furthermore, we find dynamically stable ground state densities 
with ripple-like oscillations along the radial direction. These structured ground states exist in the 
vicinity of a dynamical radial roton-like instability. 

PACS numbers: 



I. INTRODUCTION 

Dipolar Bose Einstein condensates (BECs) such as 
magnetic Cr condensates are characterized by angle- 
dependent long-range interactions [ESQ- Usually, 
the magnetic dipole-dipole interactions compete with the 
isotropic short-range s-wave contact interactions, which 
dictate the behaviors of alkali atom BECs [4| as well as of 
two-component Fermi gases such as 40 K and 6 Li HQ. In 
Cr, the s-wave scattering length a s can be tuned to zero 
through the application of an external field in the vicin- 
ity of a Fano-Feshbach resonance, thus allowing for the 
realization of a pure dipolar Bose gas @, H, @ • Although 
seemingly simple, dipolar Bose gases in which the dipoles 
are aligned along a particular laboratory axis have been 
shown to exhibit a variety of intriguing and unique fea- 
tures such as so-called roton instabilities [13, E[ and 
structured density profiles [12, 0, G3 ■ Importantly, al- 
though the strength of the magnetic dipole-dipole inter- 
actions is fixed by the electronic structure of the atoms, 
it can be tuned effectively by varying the aspect ratio of 
the external confining potential [15l . This effective 
tunability of the dipolar interactions underlies a series of 
studies of dipolar gases and can be readily understood 
intuitively. Assuming the dipoles are aligned along the 
z-axis, an effectively one-dimensional confinement that 
predominantly allows for motion along the z-axis leads 
to an attraction between the dipoles and thus collapse. 
An effectively two-dimensional confinement that predom- 
inantly allows for motion in the xy-plane, in contrast, 
leads to a repulsion between the dipoles and thus stabi- 
lization. 

This paper considers a pure dipolar Bose gas, aligned 
along the z-axis, in an oblate trapping geometry with 
cylindrical symmetry characterized by the aspect ratio 
A. In addition, a repulsive Gaussian barrier with fixed 
width b and variable height A is added along the z direc- 
tion (centered at z = 0). When the chemical potential /i 



is much larger than the barrier height A, the Gaussian 
potential serves as a small perturbation. For /i/A -C 1, 
in contrast, the fairly narrow Gaussian barrier signifi- 
cantly modifies the system behavior and approximately 
splits the dipolar condensate into two "sheets" or "pan- 
cakes" . We determine the phase and stability diagram of 
the dipolar g function of the mean-field strength 

D and the barrier height A for various aspect ratios A. 
Our key results are as follows: (i) We find an instability 
island immersed in otherwise stable regions of the D ver- 
sus A phase diagram. Our analysis of the Bogoliubov-de 
Gennes excitation frequencies and corresponding eigen- 
modes indicates that these unstable islands arise due to 
a breathing mode-like instability, (ii) We find dynami- 
cally stable density profiles of Gaussian shape with su- 
perimposed ripple-like structure. For certain parameter 
combinations near the instability line, these ripple-like 
structures exist for vanishing as well as finite Gaussian 
barrier. As the mean-field strength increases, one of the 
Bogoliubov-de Gennes frequencies with radial roton type 
character becomes soft, inducing a dynamical instablility. 

Experimentally, the instability island predicted by our 
calculations can be probed straightforwardly by follow- 
ing two different types of trajectories in the D versus 
A phase diagram: (i) The first experiment probes the 
instability by adiabatically increasing the barrier height 
from A\ to A\ c or by adiabatically decreasing the bar- 
rier height from A-i to Ai c for fixed aspect ratio A, fixed 
mean-field strength D and fixed barrier width b (see tra- 
jectories 1 and 2 in Fig. []}. The experimental signature 
of the collapse would be the onset of significant atom loss 
at A\ c and A2 C , respectively, (ii) The second experiment 
by-passes the instability island by moving from A\ to A2 
and back to A\ by simultaneously varying A and D while 
keeping A and b fixed (trajectories 3 and 4 in Fig.Q]trace 
out one of many possible paths). In the case of Cr, vary- 
ing the mean- field strength D, which is determined by the 
number of particles N and the square of the dipole mo- 



2 



20 



10 



-i 1 1 1 1 1 r- 

3 



4 

A/E 



FIG. 1: (Color online) Probing the instability island. A 
solid line separates the mechanically stable region and the 
mechanically unstable island (labeled "U") in the mean-field 
strength D versus barrier height A phase diagram for A = 8 
and b = 0.2a z . Filled circles indicate (D, A) values at which 
the k — mode becomes dynamically unstable (see text for 
details). The instability island can be probed experimentally 
by following trajectories 1 and 2 that end at the edges of 
the dynamically unstable island (labeled Ai c and Ai c ), and 
by following trajectories 3 and 4 that encircle the instability 
island. 



merit strength d, implies varying N. Thus, trajectories 
3 and 4 can be mapped out by preparing dipolar BECs 
with varying N and observing long condensate lifetimes 
for all (D,A) parameter combinations along the chosen 
trajectories. 

The remainder of this paper is organized as follows. 
Section [TT] introduces the system under study and the 
mean-field formalism. Section IIIII presents our results 
and Sec. IIVI concludes. 



II. MEAN-FIELD DESCRIPTION OF DIPOLAR 
BECS 

Our description of the aligned dipolar Bose gas with 
vanishing s-wave interactions is based on the time- 
dependent Gross-Pitaevskii equation [It], [H, [H| 



ih 



2m 



V 2 + ^ xt (f) + 



d 2 (N-i) I y^|^,t)| 2 dV]*!rwi. (i) 



where m denotes the mass of the dipoles, d the strength 
of the dipole moment, and i? the angle between the rel- 
ative distance vector f — f and the z-axis. Throughout, 
we employ cylindrical coordinates f = (p, ip, z). The ex- 
ternal confining potential T4xt(p, z), 



V cxt (p,z) = ^m(u> 2 p 2 +u 2 z z 2 ) + Aexp (- 



z 

262 



(2) 



is characterized by the transverse and longitudinal an- 
gular trapping frequencies lo p and uj z , respectively, and 



the height A and width b of the repulsive Gaussian bar- 
rier. Throughout, we fix b, i.e., we set b — 0.2a z , where 
a z = y/KJ (mu! z ). The angular trapping frequencies de- 
fine the aspect ratio A, 



A = 



(3) 



Throughout, we consider oblate systems with A = 6 — 9. 
If the Gaussian barrier is sufficiently large, the system 
consists — to lowest order — of two staggered pancakes; in 
this case, we define an effective aspect ratio A e ff, which 
is obtained by analyzing the energy spectrum for large A 
and vanishing interactions, i.e., for D = (see below). To 
quantify the strength of the mean-field term in Eq. |1| , 
we define the dimensionless quantity D, 



D 



d 2 {N -1) 



where E„ 



huj p and a p 



(4) 



Note that the 



quantity D defined in Eq. (HI differs from the definition 
used in our previous work [201 ] but agrees with that used 
in Ref. 

We consider stationary and dynamical properties of the 
dipolar Bose gas. Our stationary calculations determine 
the energetically lowest-lying gas-like state of the dipo- 
lar system. We write \l/(r, t) = exp(— ifxt/H)tpo{pi z ): 
where -00 denotes the stationary ground state solution 
and n the corresponding chemical potential. Since the 
system is cylindrically symmetric, "00 is a function of p 
and z only and not of ip. We solve the resulting two- 
dimensional Schrodinger equation by propagating a given 
initial state in imaginary time. In selected cases, we ad- 
ditionally employ a basis set expansion type approach. 
We have checked that the results of these two distinctly 
different approaches, whose implementation follows that 
outlined in Ref. [HI, [2^, Hl| , agree within our numerical 
uncertainties. 

In addition to the ground state properties, we de- 
termine the Bogoliubov-de Gennes eigenfrequencies and 
eigenmodes, which provide insights into the dynamical 
stability of the system 0, H2] • To this end, we solve the 
Bogoliubov-de Gennes eigenequations, which are derived 
within linear response theory assuming that the ground 
state ipo is perturbed by 6ijj(r,t), 



6rp(f,t) = u(f) exp(— iuit) + v* (r) exp(iwt) , 



(5) 



iteratively using the Arnoldi method [23| . Here, u(r) and 
v(r) are the so-called Bogoliubov-de Gennes functions 
and uj is the Bogoliubov-de Gennes eigenfrequency. Fol- 
lowing the approach introduced in Ref. [23[ , we solve the 
Bogoliubov-de Gennes eigenequations for oj 2 and for the 
so-called Bogoliubov-de Gennes eigenmode /(f), where 
/(f) = u{f) + v(r). Because of the cylindrical symme- 
try of the system, the tp dependence can be separated 
off, i.e., /(f) = f(p, z)h(<p) with h((p) = exp(ik<p) / y/2n . 
This implies that the Bogoliubov-de Gennes excitation 
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frequencies u> can, for a given azimuthal quantum num- 
ber k, be obtained by solving a two-dimensional eigen- 
value problem. Typically, we consider excitations with 
up to k = 4. For the discussions in Sec. IIII1 it is im- 
portant to recall that a negative oj 2 , and thus a purely 
imaginary u>, signals a dynamically unstable ground state 
ipo(r)- Correspondingly, an analysis of the Bogoliubov- 
de Gennes eigenmodes f{r) near the dynamical insta- 
bility point provides insights into the decay mechanism; 
in particular, the correction to the time-dependent den- 
sity is, to lowest order in the perturbation Sip, given by 
2cos(urf)Vo(r)/(r). 

To elucidate the energy and length scales of V cyA {p, z), 
we consider the non- interacting system (i.e., we set D — 
0). In this case, the linear Schrodinger equation sep- 
arates and the solutions are readily obtained. In the 
(p, </?)-plane, we have a two-dimensional harmonic oscil- 
lator with eigenenergies (2n p + \k\ + l)E p , where n p and 
k denote the principal and azimuthal quantum numbers, 
respectively (n p = 0, 1, • • • and k = 0,±1, • • • ). In the 
^-direction, the energy spectrum depends on two param- 
eters, the axial frequency lo z and the barrier height A 
(recall, the barrier width b is kept constant throughout 
this paper). Figure [5] shows the lowest six eigenener- 
gies, obtained by diagonalizing the Hamiltonian matrix 
expressed in terms of harmonic oscillator functions, as a 
function of A for b = 0.2a z . For A = 0, the eigenener- 
gies follow the harmonic oscillator pattern (n z + 1/2)E Z , 
where n z = 0, 1, • • • and E z — hw z . For large A, in con- 
trast, the eigenenergies appear in roughly equally spaced 
pairs [Fig. 0(a) shows three of these pairs]. We use the 
energy spacing between the lowest two pairs to define 
an effective angular frequency w Zje ff that characterizes 
the left and the right well of the double well poten- 
tial. Extrapolating to the large barrier height limit, we 
find w ZiG ff ~ 2uj z . Correspondingly, we define an effec- 
tive harmonic oscillator length a z . c s through a z ^ g = 
^Jhj(muj z ^\ leading to a ZjC ff ~ a z j\[2. While this 
analysis is approximate — our definition of a 2iG fj, e.g., as- 
sumes a harmonic potential around the minimum of the 
left well and of the right well — , it provides insights into 
how the energy and length scales associated with the z- 
degree of freedom vary with increasing barrier height. 

Lastly, we point out that the Gaussian barrier leads 
to an energetically low-lying tunneling splitting mode for 
comparatively large A and D = 0. For A = 16E Z , e.g., 
we find that the density at z — is about 25 times smaller 
than the peak density, suggesting that the left well and 
the right well are, to a good approximation, decoupled. 
Correspondingly, as shown in Fig. HJb), the lowest fre- 
quency decreases from \lo z for A = to about 0.044w z 
for A — 16E Z . This implies that the tunneling splitting 
mode for D — is, in the large A limit and for A = 6 — 9, 
energetically lower-lying than the lowest mode along the 
p-direction which has a frequency of uj p . Our discussion of 
the excitation spectrum for vanishing D presented here 
serves as a guide to understanding the Bogoliubov-de 
Gennes excitation spectra for finite D (see Sec. Mil) . 




FIG. 2: (Color online) Energetics of the one-dimensional 
double-well potential. Solid lines show (a) the eigenenergies 
and (b) the lowest eigenfrequency of the one-dimensional dou- 
ble well potential as a function of the barrier height A. The 
spectrum consists of nearly equally spaced energy levels for 
small A and nearly equally spaced energy level pairs for large 
A. For comparison, symbols in panel (b) show the eigenfre- 
quencies obtained by solving the Bogoliubov-de Gennes equa- 
tions for D — 0; excellent agreement between the two sets of 
frequencies is found. 



III. RESULTS 

This section presents the energetics, the density pro- 
files, and selected Bogoliubov-de Gennes excitation spec- 
tra and eigenmodes for trapped dipolar Bose gases. In 
particular, we detail calculations that predict an instabil- 
ity island (see Fig.[T]) and density profiles with ripple-like 
oscillations. 

To understand how the instability island for A = 8 
emerges, Figs. [Ua)-(c) show the D versus A phase di- 
agram for fixed b and three different A values, i.e., for 
A = 6, 8 and 9. Compared to Fig. [TJ the phase dia- 
gram for A = 8 [see Fig. [3](b)] shows an extended (D, A) 
parameter regime. As in Fig. [TJ mechanically stable re- 
gions of the phase diagram are separated by solid lines 
from mechanically unstable regions of the phase diagram. 
Dynamically stable regions are separated from dynami- 
cally unstable regions by symbols; filled circles, squares 
and diamonds indicate that the dynamical instability is 
triggered by a k = 0, 2 and 3 mode, respectively. The 
mechanical and dynamical instability points that sepa- 
rate the "upper" and "lower" regions of the phase dia- 
gram differ by up to about 30 %. Along these instabil- 



4 




FIG. 3: (Color online) D versus A phase diagram for fixed 
b, b = 0.2a z , and three different A values: (a) A = 6, (b) 
A = 8, and (c) A = 9. Solid lines separate mechanically stable 
regions (labeled "So" and "S>o") from mechanically unstable 
regions (labeled "U"). The instability island centered around 
A » 3.5E Z and D w 13 for A = 8 [see panel (b)] is shown on 
an enlarged scale in Fig. [1] The dynamically stable regions 
are separated from the dynamically unstable regions by sym- 
bols: Filled circles, squares and diamonds indicate that the 
dynamical instability is triggered by a k = 0, 2 and 3 mode, 
respectively. Dashed lines separate stationary ground state 
densities with maximum at p — (labeled "So") from those 
with maximum at p > (labeled "S>o"). In certain regions of 
the phase diagram, the ground state densities possess ripple- 
like oscillations (the corresponding phase boundaries are not 
indicated here; see text and Fig. I10[) . 



ity lines, the decay is triggered by modes with vanishing 
or non-vanishing azimuthal quantum number k and has 
previously been discussed for A — in Ref. [12] and for 
A = 12E Z in Ref. [2(|. In the vicinity of the instability 
island, which is centered at A w 3.5E Z and D w 13 [see 
elliptically-shaped solid lines in Fig.[3jb) and Fig. [I], the 
k = Bogoliubov-de Gennes mode becomes soft first; 
in this region, the dynamically unstable region is only 
slightly larger than the mechanically unstable region. 

Inspection of Figs. Ela)-(c) shows that the instability 
island for A = 8 emerges from the break-up of the un- 
stable region of the phase diagram for A = 6 into two 
pieces [see Fig. |3Ja)]. For A = 7 (not shown), we find 
a somewhat larger mechanically unstable island than for 
A = 8, centered around roughly the same (D, A) values 




FIG. 4: (Color online) Energetics across the instability is- 
land. Solid and dashed lines show the total energy per par- 
ticle Etot/N and the chemical potential p, respectively, for 
A = 8 and b — 0.2a z . In panel (a), the energetics are shown 
as a function of D for A = 3.5E Z ; in panel (b), the energetics 
are shown as a function of A for D — 13. 



as for A = 8. Thus, as A increases the instability is- 
land shrinks, and it is absent for A = 9. Figure [3£b) 
shows that the instability island is surrounded by a dy- 
namically stable region, suggesting that the experimental 
sequence discussed in the introduction in the context of 
Fig. [T] is feasible with present-day technology. Specifi- 
cally, for a Cr condensate with vanishing s-wavc scatter- 
ing length, confined by a trap with lu p = 2ir x 100 Hz 
and uj z = 2ir x 800 Hz, point A\ in Fig. [T] corresponds to 
N w 7200 atoms, while the dynamical instability points 
for A — 3.5E Z in Fig. 0(b) correspond to (from bottom 
to top) N « 5900, 9300 and 14200. 

In determining the D versus A phase diagrams shown 
in Figs. [3|a)-(c), our ground state calculations employed 
a grid with widths 4 and 1E Z in D and A, respectively. 
In a second set of calculations, a finer grid in D was 
employed to determine the mechanical instability point 
more accurately and to analyze the dynamical stability. 
We cannot rule out the existence of additional instability 
islands with smaller widths than our numerical resolu- 
tion. 

Figure Q] shows the energy per particle E to t /N (solid 
line) and the chemical potential fj, (dashed line) for 
A = 8 and b = 0.2a z : Fig. Ufa) shows the energet- 
ics for A = 3.5E Z (corresponding to 28E p ) as a func- 
tion of D [i.e., along a vertical cut in Fig. [21b)], while 
Fig. 0Jb) shows the energetics for D = 13 as a function 
of A [i.e., along a horizontal cut in Fig. E^b)]. For fixed 
D [see Fig. EJb)], the energy and chemical potential in- 
crease monotonically till the mechanical instability point 
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FIG. 5: Examplary integrated density profiles n(p) as a 
function of p for A = 9, b = 0.2a z and two different (D, A) 
combinations: (a) (D, A) = (40, 0) (So type density) and (b) 
(D, A) = (36, 4E Z ) (S>o type density). Both density profiles 
correspond to dynamically stable dipolar Bose gases near the 
instability line [see Fig. E2c)]. 



is reached at A » 2.85E Z . At this point, the chemical 
potential p is about two times smaller than the barrier 
height. Beyond the mechanically unstable region, i.e., 
for A > 4AE Z , both E tot /N and p continue to increase 
monotonically. The energy per particle and the chemical 
potential behave similarly as a function of D for fixed A 
[see Fig.UJa)]. Interestingly, the presence of the instabil- 
ity island manifests itself in the form of a "gap" in the 
energy and the chemical potential but does not notably 
change the slope of either of these quantities (i.e., E tot /N 
and p appear to be changing smoothly across the insta- 
bility island). We find a similar behavior for a subset of 
the Bogoliubov-de Gennes eigenmodes (see below). 

In addition to the energetics, we analyze the structural 
properties of dipolar Bose gases in a double-well poten- 
tial. Dashed lines in Figs.[3]Ja)-(c) separate the parameter 
region where the ground state density IV'ol 2 has its max- 
imum at p — (labeled So) from that where [f/'ol 2 has its 
maximum at p > (labeled S>o)- The transition from 
the So type to the S>o type density is smooth and the 
dashed lines are based on an analysis of the integrated 
density n(p), n(p) — 2ir f | -i/? (/°, z)\ 2 dz. We denote the p 
value at which n(p) is maximal by p max and define the 
transition from So to S>o type densities, i.e., the dashed 
lines in Fig. [H by the condition n(0)/n(p max ) = 0.98. 
Figures [5Ja) and (b) show examplary integrated densi- 
ties n(p) of types So and S>o, respectively. The latter 
density profile has, for a system with vanishing barrier, 
been previously termed "red blood cell" jl2j; this name 
is motivated by the fact that the system's isodensity plot 
resembles the shape of a red blood cell. For all parame- 
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FIG. 6: (Color online) Bogoliubov-de Gennes eigenspec- 
trum as a function of D for k = 0, A = 3.5E Z , X = 8 and 
b — 0.2a z . Solid lines show frequencies for which the corre- 
sponding eigenmodes / have nodal lines that are parameter- 
ized by p w constant. Dashed lines show frequencies for which 
the corresponding eigenmodes have, in addition to other nodal 
lines, a nodal line given by z — 0. Filled circles and squares 
mark those D values for which Fig. [7] and Fig. [8] respectively, 
show eigenmodes. 



ter combinations investigated, the dynamical instability 
is triggered by a k = mode if the ground state density 
is of So type and by a finite k mode if the ground state 
density is of S>o type. In particular, Fig. 0(b) shows 
that the density profiles in the vicinity of the instability 
island have a simple Gaussian shape and that the dynam- 
ical instability in this regime is associated with a k = 
mode. 

To shed light on the dynamics in the vicinity of the in- 
stability island, Fig. [5] shows the Bogoliubov-de Gennes 
excitation spectrum as a function of D for k = 0, 
A = 3.5E Z , A = 8 and b — 0.2a z . The lowest excitation 
frequency becomes purely imaginary at D sa 10, signal- 
ing the dynamical instability. The frequency "recovers" 
at somewhat larger D values (D ss 16). We find that the 
real parts of the eigenfrcquencies with non- vanishing az- 
imuthal quantum number k remain finite in the vicinity 
of the instability island, implying that the dynamics in 
this parameter space of the phase diagram is dominated 
by k — modes. Our analysis of the Bogoliubov-de 
Gennes eigenmodes (see Figs. [7] and |S] as well as the dis- 
cussion below) shows that the eigenmodes / correspond- 
ing to the Bogoliubov-de Gennes excitation frequencies 
shown by solid lines in Fig. [6] have nodal lines that are to 
a good approximation independent of z while those cor- 
responding to the Bogoliubov-de Gennes excitation fre- 
quencies shown by dashed lines have, among other nodal 
lines, a nodal line given by z = 0. The eigenfrequen- 
cies shown by solid lines are affected by the instability 
(i.e., the slope of these frequencies changes near the dy- 
namical instability), while those shown by dashed lines 
appear to change smoothly across the instability. This 
suggests that these two classes of eigenfrequencies are 
approximately decoupled. We find that the finite k fre- 
quencies also change as though they are uneffected by 
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FIG. 7: (Color online) Examplary Bogoliubov-de Gennes 
eigenmodes with breathing mode character for k — 0, A — 
3.5E Z , A = 8, b = 0.2a z and (a) D = 6, (b) D = 9.75, (c) 
D — 17, and (d) D — 20. The corresponding Bogoliubov-de 
Gennes eigenfrequencies are shown by circles in Fig. [6] and in 
Fig. [5] The contours are chosen equidistant, with solid and 
dashed lines corresponding to positive and negative values of 
/. The dash-dotted lines indicate the nodal lines of /. 
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FIG. 8: (Color online) Examplary Bogoliubov-de Gennes 
eigenmodes f(p, z) with tunneling splitting character for k = 
0, A = 3.5£ z , A = 8, b = 0.2a z and (a) D = 9.75 and (b) 
D — 17. The corresponding Bogoliubov-de Gennes eigenfre- 
quencies are shown by squares in Fig. [6] The contours are 
chosen equidistant, with solid and dashed lines corresponding 
to positive and negative values of /. The dash-dotted lines 
indicate the nodal lines of /. 



the instability island. For larger D values (D > 24.5), 
however, modes with finite azimuthal quantum number 
become relevant in determining the system's dynamical 
stability for a fairly wide range of A values. In fact, as 
indicated in Fig. [3Jb), the instability for comparatively 
large D values and A = to A w 10E Z is triggered by 
k = 2 and 3 modes (i.e., for these A values the energet- 
ically lowest-lying k = 2 or 3 frequency becomes purely 
imaginary before the lowest k = frequency does). 

Figure[5]shows the real part of the energetically lowest- 
lying k = Bogoliubov-de Gennes eigenfrequency as a 
function of D for A = 8, b — 0.2a z and various A val- 
ues. The frequency shown in Fig.[SJc) for A = 3. 5E Z has 
just been discussed in the context of Fig. [Sj Inspection 




FIG. 9: (Color online) Real part of the energetically lowest- 
lying Bogoliubov-de Gennes eigenfrequency lu as a function 
of D for k — 0, A = 8, b = 0.2a z and five different A values: 
(a) A = 0, (b) A = 2E Z , (c) A = 3.5E Z , (d) A = 5E Z , and (e) 
A — 7E Z . Filled circles in panel (c) mark those D values for 
which Fig. [7] shows eigenmodes. 



of Fig. [5] shows that the vanishing of the real part of 
the energetically lowest-lying breathing mode type fre- 
quency for A = 3.5E Z is accompanied by a "dip" in the 
corresponding frequencies for smaller and larger A [see 
Figs. [UJb) and (d)] around D ss 10 and 15, respectively. 
For these A values, however, the mode "recovers" as D 
increases and the system is mechanically and dynamically 
stable up to comparatively large D values. 

We now discuss the Bogoliubov-de Gennes eigenmodes 
shown in Fig. [7] and in Fig. [8] in more detail. The k = 
eigenmodes shown in Fig. [H] show one nodal line that is 
given by z = and a second nodal line that depends on 
p and z. The former nodal line reflects the symmetry 
of the confining geometry and corresponds to oscillations 
between the left well and the right well, while the lat- 
ter is of secondary importance since the amplitude of / 
is small along this nodal line. We classify the eigenfre- 
quencies associated with eigenmodes of the type shown 
in Fig. [5] as having tunneling splitting mode character. 
The k = eigenmodes shown in Fig. in contrast, are 
characterized by a single nodal line which depends to first 
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order only on p and not on z. We classify the eigenfre- 
quencies associated with eigenmodes of the type shown 
in Fig. [7| as having breathing mode character. An analy- 
sis of the eigenmodes corresponding to the energetically 
lowest-lying k — frequency for other D values but the 
same A, A and b shows that the nodal line is located 
at p = a p for D = and bends slightly as D increases 
while remaining located around p m a e [see Fig.[7]Ja)-(b)]. 
Past the instability, the nodal line of / is again located at 
p » a p [see Fig. EJc)] and moves to larger p values with 
increasing D [see Fig. EJd)] . Eventually, a new nodal line 
moves in from p = 0. Near the dynamical instability at 
D w 36, the eigenmode is characterized by two approx- 
imately equally spaced nodal lines with approximately 
constant p (not shown; qualitatively, the eigenmode is 
similar to those shown in Fig. [TT] but without the third 
nodal line). For large D [D sa 36, see Fig. EJc)], the 
instability thus acquires some similarities with a radial 
roton instability (see below). 

Lastly, we show that cylindrically symmetric systems 
with modest A support, in addition to density profiles 
of simple Gaussian shape and of red blood cell type 
shape, density profiles with ripple- like oscillations. Fig- 
ures Hnia)-(c) show examplary integrated density profiles 
n(p) with ripple- like oscillations for A = 9 and three dif- 
ferent (D,A) parameter combinations. In general, we 
find density profiles with ripple-like oscillations in a rela- 
tively small region near the dynamical instability points 
for A = 9 but not for A = 8 and 6 [the parameter combi- 
nations that support ripple-like density profiles fall into 
the So region of the phase diagram shown in Fig. [3Jc); 
they are not indicated explicitly in Fig. EDJc)]. Ripple- 
like structures have very recently been predicted to ex- 
ist in dipolar systems with non-vanishing s-wave scat- 
tering length consisting of pancake-shaped stacks, where 
each pancake is characterized by a large aspect ratio A, 
e.g., A = 340 [H|. Ripple- like oscillations have also been 
predicted to exist for dipolar systems with vortices [l4T |. 
Here, we find ripple-like oscillations in the dynamically 
stable region of the phase diagram of pure dipolar gases 
confined by a cylindrically-symmetric trapping geometry 
with and without Gaussian barrier and with modest as- 
pect ratios A. 

The integrated density profiles shown in Fig. fTUT aVfc) 
for A — to 12E Z possess density oscillations with char- 
acteristic length scale of the order of a p . Figures fTTTa)- 
(c) show the corresponding eigenmodes / of the energet- 
ically lowest-lying k = Bogoliubov-de Gennes eigen- 
frequency. The Gaussian barrier modifies the amplitude 
of the eigenmode near z = but otherwise leaves the 
overall structure of the eigenmodes unaffected. All three 
eigenmodes possess roughly equally spaced nodal lines 
that are to a good approximation independent of z, and 
which can, roughly speaking, be characterized — just as 
the ripple- like oscillations — by p n odo ~ na p , where n = 1, 
2 and 3. Using A = 9, it can be seen that Pnode is deter- 
mined by p nodc w 7ra p /\/A or p nodc w na z pH E2- The 
latter expressions reflect that the condensate develops 




FIG. 10: Integrated ground state density profiles n(p) with 
ripple-like oscillations as a function of p for A = 9, b = 0.2a z 
and different (D, A) combinations: (a) (D, A) — (50,0), (b) 
(D,A) = (68,8£ z ) and (c) (D, A) = (99, 12E Z ). The length 
scale at which the ripple-like oscillations occur is to a good 
approximation independent of the barrier height A (see text 
for details). All three density profiles correspond to dynam- 
ically stable dipolar Bose gases near the instability line [see 
Fig. EJc)]. The eigenmodes / of the corresponding energet- 
ically lowest-lying k — Bogoliubov-de Gennes frequencies 
are shown in Fig. 1111 



three-dimensional character, i.e., that it "gets chopped 
up" into smaller pieces along the p-dircction. The dy- 
namical instability that arises when D is increased some- 
what compared to the D values chosen in Figs. ITD1 and [TT1 
[see also Fig. [He)] is thus identified as a radial roton-likc 
instability [lCJ, ll2| . Intuitively, one might expect that the 
length scale that characterizes the ripple-like oscillations 
and the nodal lines of the Bogoliubov-de Gennes eigen- 
modes would depend on the barrier height A or, equiva- 
lently, the effective aspect ratio A e fj. Figures [TU1 and [TTl 
however, suggest that the system behavior is determined 
by A instead. We note, though, that the change of A e ff 
as a function of A is fairly small, implying that further 
studies are needed to determine the relevant length scale 
unambiguously. 



IV. SUMMARY 

We have investigated oblate dipolar Bose gases with 
vanishing s-wave scattering length in a double well po- 
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FIG. 11: (Color online) Examplary Bogoliubov-de Gennes 
eigenmodes f(p,z) for k = 0, A = 9, b = 0.2a z and differ- 
ent (D,A) combinations: (a) (D, A) = (50,0), (b) (D,A) = 
(68,8£ z ) and (c) (D, A) = (99, 12E Z ). The contours are cho- 
sen equidistant, with solid and dashed lines corresponding 
to positive and negative values of /. The dash-dotted lines 
indicate the nodal lines of /. The corresponding integrated 
ground state densities n(p) are shown in Fig. [TUJ 

tential as a function of the barrier height A, the mean- 
field strength D and the aspect ratio A. Our stationary 
and dynamical mean-field calculations add two new as- 
pects to the already extensive list of intriguing behaviors 
of dipolar Bose gases: (i) We find an instability island im- 
mersed in a mechanically and dynamically stable region 
of the phase diagram. This instability island can, as has 
been outlined in Sees. Uland lllli be probed with present- 
day technology and arises due to the going soft of a radial 
breathing mode-like frequency, (ii) We find structured 
ground state densities with ripple-like oscillations in the 



dynamically stable region of the phase diagram for mod- 
erate aspect ratios and cylindrically-symmetric confining 
geometries. These ripple-like oscillations exist, for ap- 
propriately chosen parameter combinations, for vanishing 
and non-vanishing Gaussian barriers; an increase of the 
mean-field strength eventually induces a radial roton-likc 
instability. In the future, it will be interesting to inves- 
tigate if the ground state densities with ripple-like oscil- 
lations can be probed experimentally using time-of-fiight 
expansion techniques. Theoretically, this question can 
be addressed by determining the real time dynamics af- 
ter turning off the confining potential and by comparing 
the results with those for simple Gaussian density pro- 
files. We emphasize that the emergence of the instability 
island and of density profiles with ripple-like oscillations 
is unique to dipolar gases — these features are a direct 
signature of the long-range, anisotropic dipole-diple in- 
teraction and absent for s-wave interacting Bose gases. 

Current experiments are directed at loading molecular 
gases with appreciable electric dipole moment [25|, (2(| 
into one-dimensional lattices. The premise of these ex- 
periments is to exploit the effectively two-dimensional ge- 
ometry of each lattice site as a stabilization mechanism. 
In this context, several theoretical studies 0, H3, [H[ 
have investigated stacks of effectively two-dimensional 
dipolar Bose gases. Common to these studies is the as- 
sumption of a very tight confinement along the z-axis, 
i.e., a large aspect ratio for each individual lattice site, 
and a comparatively large spacing between individual lat- 
tice sites. In those studies, roton and phonon instabili- 
ties emerge. The studies presented here can be viewed 
as an alternative approach to tackling the multi-site op- 
tical lattice system. Using a single mean-field equation 
and working in a regime where the extent of the cloud 
in one well is approximately equal to the separation be- 
tween clouds, the full system dynamics is investigated 
for moderate aspect ratios. Future studies will investi- 
gate the system behavior as a function of the separation 
between sites, using a single mean-field equation and rel- 
atively small barrier heights as well as a set of coupled 
mean-field equations; these future studies will facilitate 
direct comparisons between the approach pursued here 
and Refs. [mEHH- 
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